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Parallel  channels  have  many  advantages,  such  as  low  pressure  drop  and  easy  fabrication,  but  they  may 
cause  flow  maldistribution  which  would  result  in  low  reaction  efficiency.  This  study  presents  an  analytical 
model  to  calculate  the  flow  distribution  of  the  parallel  channels  based  on  the  assumption  of  the  analogy 
between  fluid  flow  and  electrical  network.  The  model,  which  ultimately  releases  from  the  solution  of  a  set 
of  nonlinear  equations,  is  validated  by  comparing  with  the  results  obtained  from  three-dimensional  com¬ 
putational  fluid  dynamics  (CFD)  simulations.  Consequently,  the  model  is  used  to  optimize  the  geometric 
dimension  of  a  parallel  plate  to  obtain  a  uniform  flow  field  distribution. 

©  2009  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

Proton  exchange  membrane  fuel  cells  (PEMFC)  are  considered  to 
be  advanced  power  sources  for  their  high  efficiency  and  low  pollu¬ 
tion.  Bipolar  plates  are  important  components  of  the  PEMFC,  which 
help  to  distribute  the  reactant  gas,  collect  current,  provide  structure 
support,  facilitate  water  and  heat  management  [1].  Many  designs 
of  flow  configurations,  for  example,  serpentine,  parallel  discontin¬ 
uous  [2]  and  interdigitated  channel  (Fig.  1),  are  currently  used. 
There  are  two  principal  considerations  in  the  choice  of  a  partic¬ 
ular  configuration:  (i)  the  overall  pressure  drop;  (ii)  the  degree  of 
non-uniformity  of  flow  distribution  over  the  plate  [3].  Serpentine 
channels  have  good  flow  distribution  but  may  have  high  pressure 
drop,  while  parallel  channels  have  low  pressure  drop  but  may  cause 
flow  maldistribution  which  would  result  in  severe  problems.  For 
example,  some  channels  may  be  staved  of  reactants,  while  others 
may  have  excessive  reactants. 

Generally,  there  are  two  kinds  of  method  calculating  the  flow 
distribution  of  bipolar  plates:  (i)  CFD  simulation  and  (ii)  analytical 
method.  With  the  development  of  computer  technology,  CFD  sim¬ 
ulations  are  usually  more  accurate  than  analytical  models,  but  they 
always  cost  more  time  and  money.  As  a  result,  it  is  hard  to  opti¬ 
mize  the  geometric  parameters  of  the  bipolar  plates  to  fix  the  flow 
maldistribution  of  parallel  channels  using  CFD  simulations  due  to 
numerous  calculations. 
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Some  researches,  dealing  with  the  flow  maldistribution  of  the 
parallel  channels,  have  been  reported.  Ma  et  al.  [4]  considered  the 
consequences  of  a  non-uniform  gas  distribution  using  a  computa¬ 
tional  model.  Pei  [5]  studied  the  hydrogen  pressure  characteristics. 
Kee  et  al.  [6]  raised  the  possibility  of  maldistribution  and  presented 
a  numerical  model  to  calculate  the  flow  distribution  in  a  Z-type 
parallel  configuration.  Maharudrayya  et  al.  [7]  analyzed  U-type 
parallel  configurations  and  obtained  closed  form  analytical  solu¬ 
tions  to  calculate  the  flow  distribution  and  the  pressure  drop  in 
simple  parallel  channels.  Jung  et  al.  [8]  developed  a  2-D  model  to 
analyze  the  flow  and  pressure  distributions.  Shimpalee  et  al.  [9] 
numerically  investigated  how  serpentine  flow-fields  with  different 
channel/rib’s  cross-section  areas  affect  performance  and  species 
distributions  for  both  automotive  and  stationary  conditions.  Fur¬ 
thermore,  Koh  et  al.  [10]  presented  a  resistance  concept  based 
model  for  to  evaluate  the  flow  distribution  among  cells  in  a  fuel 
cell  stack.  Ma  [11]  and  Karimi  [12]  used  network  flow  analysis  to 
calculate  flow  distribution  within  fuel  cells  with  nonlinear  meth¬ 
ods.  However,  a  simple  design  procedure  that  can  be  used  for  the 
optimization  of  the  parallel  flow  field  is  not  available. 

Against  this  background,  the  present  work  proposes  a  simple 
linear  computational  model  for  the  purpose  of  optimization  design 
to  resolve  the  problem  of  maldistribution.  Three-dimensional  CFD 
simulations  are  also  performed  to  validate  the  model  and  optimiza¬ 
tion  results. 

2.  Theory 

2.1  Assumptions 

The  model  is  set  up  based  on  the  following  assumptions: 
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Nomenclature 

Ac  cross-sectional  area  of  channel  (mm2) 

Ah  cross-sectional  area  of  header  (mm2 ) 

At  cross-sectional  area  of  the  ith  channel  (mm2 ) 

bc  depth  of  channel  (mm) 

Dc  hydraulic  diameter  of  channel  (mm) 

/  friction  factor 

F  non-uniformity  index 

Lc  length  of  channel  (mm) 

l(i)  width  of  channel  (mm) 

L(z)  length  of  the  header  between  each  channel  (mm) 

n  number  of  channels 

Pin  pressure  at  inlet  (Pa) 

Pout  pressure  at  outlet  (Pa) 

Pc  cross-sectional  perimeter  of  channel,  mm 
r2-  resistance  of  the  ith  channel 

R  resistance 

Ri  resistance  of  the  inlet  header  between  ith  and  i  +  1th 

channel 

Rtf  resistance  of  the  outlet  header  between  ith  and 

i  +  1th  channel 

Ri  resistance  of  the  center  of  L  shape  corner 

Rj  resistance  of  the  center  of  T  shape  corner 

Re  Reynolds  number 

v  mean  velocity  (mm  s-1) 

Vi  mean  velocity  in  the  ith  channel  (mm  s-1 ) 

Vf  mean  velocity  in  the  inlet  header  between  ith  and 

i  +  1th  channel  (mms-1) 

Vif  mean  velocity  in  the  outlet  header  between  ith  and 
i  +  1th  channel  (mms-1) 

V[n  inlet  velocity  of  the  plate  (mm  s-1 ) 
wc  width  of  channel  (mm) 

Greek  letters 

a  channel  aspect  ratio  (width/depth) 

ft  viscosity 

p  mass  density 

rw  wall  shear  stress 


1.  The  analysis  here  assumes  the  flow  properties  (e.g.  the  mass 
density,  the  viscosity)  are  constant. 

2.  The  Reynolds  numbers  are  low,  the  flow  is  considered  as  laminar 
flow  and  the  temperature  is  nearly  uniform. 

3.  The  reactant  gas  is  considered  to  flow  in  impervious  channels.  It 
does  not  consider  the  mass  flow  between  the  reactant  channel 
and  electrolyte.  It  is  a  significant  assumption  which  is  necessary 
for  simplifying  the  model  [13]. 

2.2.  Hagen-Poiseuille  flow 

For  a  straight  channel,  steady  Hagen-Poiseuille  flow,  pressure 
drop  is  needed  to  balance  shear  stress  caused  by  the  walls. 

(Pin  —  Pout  Me  =  t\N PcLc  (1) 


where  v  is  the  mean  velocity  of  the  channel.  Here, /is  the  friction 
factor  given  by  the  empirical  correlation  of  Kays  and  Crawford  [14] 

Re/  =  13.84  +  10.38 exp  (~~ )  (3) 


where  a  =  wc/bc  is  the  channel  aspect  ratio.and  the  Reynolds  num¬ 
ber  is 


(4) 


where  Dc  is  the  hydraulic  diameter  Dc  =  4wcjbc/2(wc  +  bc). 
From  Eqs.  (1  )-(4),  we  get 


A?  =  Pin  -  Pout  = 


1  (Re/^PcIc-  p- 

2  DCAC  V  =  RV 


where 


1  (Re/)/xPcLc 

2  DcAc 


(5) 

(6) 


R  can  be  considered  as  the  fluidic  resistance  of  the  channel.  Eq. 
(5)  shows  that  for  a  given  dimensions  straight  channel,  the  pres¬ 
sure  drop  is  proportional  to  the  mean  velocity  because  the  fluidic 
resistance  of  channel  is  constant. 


2.3.  Resistance  at  corners 

The  flow  resistance  at  the  straight  channels  is  discussed  above, 
and  the  resistance  at  the  corners  will  be  discussed  in  this  section.  As 
shown  in  Fig.  2,  there  are  two  kinds  of  corners,  T  shape  and  L  shape, 
for  a  parallel  channel.  Rit  Az,  P2  are  the  flow  resistance,  the  cross- 
sectional  area  and  the  cross-sectional  perimeter  of  each  straight 
channel,  respectively.  R2  can  be  calculated  by  Eq.  (6)  if  geometric 
dimensions  are  given.  The  arrangement  strategy  of  the  resistance 
at  these  places  will  be  discussed  below. 

2.3.1.  T  shape 

The  discrete  model  for  T  shape  corner  is  shown  in  Fig.  2(c),  RT 3  is 
quite  small  compared  with  R3  because  of  the  long  length  of  channel. 
Then  RT3  can  be  neglected. 

Rj 2  can  be  considered  according  to  the  force  balance: 


PiAh  -  p2Ah  =  x wAr 


(7) 


where  AT  represents  the  contact  surface  at  T  shape  corner.  For  a  T 
shape  corner,  Aj  is  shown  in  Fig.  3  with  black  color. 

Hence,  RT2  is  obtained  with  Eqs.  (2)-(4)  and  (7) 


Rj2 


M(R  e/Mx 
2  DcAh 


(8) 


2.3.2.  L  shape 

The  discrete  model  for  L  shape  corner  is  shown  in  Fig.  2(D),  Ri  is 
quite  small  compared  with  R2  because  of  the  long  length  of  channel. 
Then  RL  can  also  be  neglected. 

After  the  rearrangement,  the  resistance  at  the  corners  will  not  be 
involved  into  the  resistance  matrix,  because  it  has  been  considered 
in  the  calculation  of  straight  channels.  This  rearrangement  facili¬ 
tates  the  calculation  and  the  availability  of  the  rearrangement  will 
be  validated  in  Section  3. 


where  pin  and  pout  are  the  pressure  of  the  inlet  and  outlet  of  the 
channel,  respectively.  Ac  is  the  channel  cross-sectional  area.  Pc  is 
the  channel  cross-sectional  perimeter,  and  Lc  is  the  channel  length. 

The  shear  stress  can  be  represented  by  using  the  friction  factor 
/as  follows: 


f  = 


tw 

1/2  pp- 


(2) 


2.4.  Analytical  model  for  parallel  channel 

Generally  speaking,  there  are  two  types  of  parallel-channel  con¬ 
figurations,  namely  Z-type  and  U-type  (Fig.  4).  Each  type  has  an  inlet 
header  and  an  outlet  header.  In  a  Z-type  configuration,  the  inlet  is 
near  the  first  channel  and  the  outlet  is  near  the  last  channel.  While 
in  a  U-type  configuration,  both  the  inlet  and  outlet  are  near  the  first 
channel. 
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Fig.  2.  (A)  T  shape  comer,  (B)  L  shape  corner  (C)  discrete  model  for  T  shape  corner 
and  (D)  discrete  model  for  L  shape. 


In  order  to  get  velocity  distribution  of  a  parallel  configuration 
with  n  channels,  a  parallel  configuration  with  3  channels  will  be 
discussed  first. 

2. 4 A.  Model  for  a  3-channel  plate 

2 AAA.  Z-type  configuration.  Consider  a  Z-type  parallel  configura¬ 
tion  with  3  channels  shown  in  Fig.  5A.  The  flow  field  is  discretized 
first  as  shown  in  Fig.  5B,  and  the  analytical  model  will  be  set  up 


Fig.  3.  Schematic  diagram  of  AT. 


based  on  solving  the  pressure  and  mass  balance  equations.  It  is 
assumed  that  the  cross-sectional  area  of  the  inlet  header,  which 
is  represented  by  Ah,  is  uniform,  and  is  equal  to  that  of  the  outlet 
header,  the  cross-sectional  area  of  the  ith  channel  represented  by 
Ai  is  not  necessary  uniform. 

As  shown  in  Fig.  5B,  Pz  is  the  pressure  at  dividing  or  combing 
points,  Vj  is  the  mean  velocity  in  the  inlet  header  between  ith 
and  i  +  lth  channel,  Vy  is  the  mean  velocity  in  the  outlet  header 
between  ith  and  i  +  lth  channel,  ft*  is  the  resistance  of  the  inlet 
header  between  ith  and  i  +  lth  channel,  ffy  is  the  resistance  of  the 
outlet  header  between  ith  and  i  +  lth  channel,  and  rz  are  the  mean 
velocity  and  the  resistance  of  the  ith  channel,  respectively. 

According  to  the  pressure  balance, 


v\  n  +vvrv 

=  Vifti  + 1^2  Tz 

O) 

V2r2  +  V'l'R'l' 

=  V2R2  +  P3r3 

(10) 

Inlet 


Outlet  °utlet 


Fig.  4.  Schematic  diagram  of  (A)  Z-type  and  (B)  U-type  parallel-channel  flow  configurations. 
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Inlet  header 


Fig.  5.  Schematic  diagram  of  (A)  Z-type  configuration  and  (B)  discrete  channel  model. 


According  to  the  mass  balance, 


Mh  =  Mi+Vi4  (li) 

V1Ah  =  V2Ah  +  v2A2  (12) 

V2Ah  =  Ms  (13) 

Vi+Vy^Vin  (14) 

V2  +  V2'=Vin  (15) 

Combining  Eqs.  (9)-(15), 

M(3)  =  R(3)V(3)  (16) 

where 

~-VmRv~ 

-VinR2' 


2.4A.2.  U-type  configuration.  Similarly,  consider  a  U-type  parallel 
configuration  with  3  channels  as  shown  in  Fig.  6. 

According  to  the  pressure  balance, 


V\  n  =  V\R\  +v2r2  +  VvRv 

(17) 

V2r2  —  V2r2  ~\~  P3P3  +  ^2,R2' 

(18) 

According  to  the  mass  balance, 

^inAi  =  V\A\  +  V\Ah 

(19) 

V\  Ah  =  V2Ah  +  v2A2 

(20) 

V2Ah  =  v3A3 

(21) 

V 1  =  Vv 

(22) 

£ 

11 

£ 

(23) 

Combining  Eqs.  (17)-(23), 

M(3)  =  R(3)V(3) 

where 


Tl 

-r2 

0 

-(ri+rv) 

0 

-  0  - 

0 

r2 

-r3 

0 

~(R2  +R2') 

0 

R0)  = 

Ai 

0 

0 

Ah 

0 

(16b) 

M(3)  = 

VinAh 

0 

a2 

0 

-Ah 

Ah 

0 

_  0 

0 

a3 

0 

-Ah 

0 

~V\ " 

-r2 

0 

~{R\+Rv) 

0 

v2 

0 

*2 

-r3 

0 

~{R2  +R2') 

VC3)  = 

v3 

(16c) 

R(3)  = 

Ai 

0 

0 

Ah 

0 

v 1 

0 

A2 

0 

-Ah 

Ah 

y2. 

.  0 

0 

A3 

0 

-Ah 

M(3)  = 


VfcA 

0 

0 


(24) 


(24a) 


(24b) 


Inlet  header 


Fig.  6.  Schematic  diagram  of  (A)  U-type  configuration  and  (B)  discrete  channel  model. 
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V(3)  = 


V\ 

V2 

v3 

V 1 

LV2 


(24c) 


Eqs.  (16)  and  (24)  inform  that,  for  a  parallel  configuration  with 
3  channels,  rit  ffy  and  Rf  is  known  according  to  Eqs.  (6)  and  (8).  If  we 
know  the  inlet  velocity  Vin,  then,  we  can  solve  the  matrix  equation 
to  get  the  mean  velocity  of  each  channel. 


V(n)  = 


V\ 

V2 

Vn 

V 1 
V2 

Vn- 1 


(25c) 


2.4.2.  Model  for  an  n-channel  plate 

Similarly,  we  now  consider  a  parallel  configuration  with  n  chan¬ 
nels. 


2.4.2 A.  Z-type  configuration.  The  schematic  diagram  of  the  contin¬ 
uous  channel  and  discrete  channel  model  of  Z-type  configuration 
are  shown  in  Fig.  7. 

Like  the  method  mentioned  in  Section  2.4.1,  a  matrix  equation 
similar  to  Eq.  (16)  will  be  set  up. 


M(n)  =  R(n)V(n) 


(25) 


The  dimension  of  resistance  matrix  R  is  expanded  from  5  x  5  to 
(2n  -  1)  x  (2n  -  1). 


2.4.22.  U-type  configuration.  The  schematic  diagram  of  the  contin¬ 
uous  channel  and  discrete  channel  model  of  U-type  configuration 
are  shown  in  Fig.  8. 

A  matrix  equation  for  U-type  is  set  up  by  the  method  mentioned 
above. 


M(n)  =  R(n)V(n) 

where 


M(n)  = 


'  — Vin  RV  ' 

"  0  > 

— VinR2 , 

0 

>  n  -  1 

^in^(n-iy 

0, 

Vin^h 

(25a) 

M(n)  = 

0  > 

Mn^h 

0  > 

0 

>n-  1 

0 

>n-  1 

_  0, 

_  0  , 

(26) 


(26a) 


where 


r\  ‘ 

-;2 

0 

0  0  ...O' 

O 

+ 
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...o 

r\ 

~r2 

0 

0  0  ...0 
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+ 

1 
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o 
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^2 

~r3 

0  0  ...0 
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-(4+4) 

0  0 
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(26b) 
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Inlet  header 


(A) 


Outlet  header 


Fig.  7.  Schematic  diagram  of  (A)  Z-type  configuration  and  (B)  discrete  channel  model. 


Inlet  header 


Fig.  8.  Schematic  diagram  of  (A)  U-type  configuration  and  (B)  discrete  channel  model. 
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pi 


3.  Validation 


V(n)  = 


V2 

Vn 
V i 
V2 

Vn- 1 


(26c) 


The  dimension  of  resistance  matrix  ft  is  also  expanded  from  5x5 
to  (2n  -  1)  x  (2n  - 1). 

As  shown  in  Eqs.  (25)  and  (26),  there  are  2n  -  1  unknown  vari¬ 
ables  and  2n  - 1  equations.  So  we  can  solve  the  matrix  equations  to 
get  the  velocity  of  each  channel. 


The  accuracy  of  the  presented  model  is  verified  by  comparing 
the  result  of  the  discrete  model  with  that  of  the  three-dimensional 
CFD  simulations. 

Here,  the  flow  distribution  of  parallel  channel  configurations  are 
simulated  using  the  commercial  CFD  code  CFD-ACE+  developed  by 
CFD  Research  Corporation. 

First,  we  study  the  flow  distribution  of  11 -channel  plates.  Then, 
the  number  of  channels  is  extended  to  21  to  examine  the  accu¬ 
racy  of  the  analytical  model  with  more  channels.  The  geometrical 
dimensions  of  the  plate  are  shown  in  Fig.  9. 

As  shown  in  Fig.  9,  channels  and  ribs  are  distributed  uniformly 
with  a  width  of  1.5  mm,  and  the  width  of  the  inlet  header  and  outlet 
header  are  3  mm,  both  of  the  depth  of  the  channels  and  headers  are 
0.6  mm. 

A  comparison  of  the  velocity  distribution  obtained  form  the  dis¬ 
crete  model  and  CFD  simulation  is  shown  in  Fig.  10.  There  is  an 
acceptable  agreement  between  the  analytical  and  CFD  results. 


(A) 


(C) 


Fig.  9.  Schematic  diagrams  of  (A)  Z-type  11-channel  (B)  Z-type  21-channel  (C)  U-type  11-channel  and  (D)  U-type  21-channel  plate  (units:  mm). 
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Fig.  10.  Comparison  of  the  results  from  CFD  and  analytical  model  (A)  Z-type  11-channel  (B)  Z-type  21-channel  (C)  U-type  11-channel  and  (D)  U-type  21 -channel  plate. 


4.  Optimization 

As  mentioned  above,  parallel  channels  have  low  pressure  drop 
but  they  also  cause  flow  maldistribution  which  can  also  be  observed 
from  Fig.  10.  It  is  possible  to  have  a  uniform  flow  distribution  in 
the  parallel  channels  by  having  a  very  large  header  cross-sectional 
area  compared  with  that  of  the  channel.  Unfortunately,  this  will 
decrease  the  active  geometrical  area  of  the  bipolar  plate.  A  uniform 
flow  distribution  in  the  parallel  channels  can  also  be  obtained  by 
changing  the  cross-sectional  area  of  the  channels,  increasing  the 
cross-sectional  area  of  the  central  channels  which  has  low  veloc¬ 
ities,  decreasing  the  cross-sectional  area  of  the  channels  near  the 
inlet  and  outlet  which  has  high  velocities. 


According  to  the  model  proposed  in  section  3,  high-performance 
language  MATLAB  program  is  developed  to  optimize  the  flow 
distribution  of  an  11 -channel  Z-type  plate.  The  original  geomet¬ 
ric  parameters  of  the  plate  have  been  shown  in  Fig.  9(A)  and 
the  variables,  object  function  and  constraints  will  be  discussed 
below. 

4.1  Variables 

The  geometric  dimensions  of  the  plate  are  treated  as  constants 
except  the  width  of  each  channel  /(i)  and  the  length  of  the  header 
between  each  channel  L(z)  (Fig.  11).  For  an  11 -channel  plate,  there 
are  21  variables:  /( 1),  /( 2), . . .,  1(11)  and  1(1),  1(2), . . .,  1(10). 


Inlet  header 


Outlet  header 


Fig.  11.  Variables  of  the  bipolar  plates. 


Table  1 

Widths  of  the  channels  (units:  mm). 
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1(1) 

1(2) 

1(3) 

1(4) 

f(5) 

1(6) 

1(7) 

1(8) 

1(9) 

1(10) 

1(11) 

0.8 

0.8750 

0.9514 

1.0237 

1.0779 

1.0983 

1.0779 

1.0237 

0.9514 

0.8750 

0.8 

4.2.  Objective  function 

Non-uniformity  index,  F,  is  defined  as  follows  [6,7]: 

y/TTan-vf 

F  =  — - - - 

nv 


Table  2 

Lengths  of  the  headers  (units:  mm). 


1(1) 

1(2) 

1(3) 

1(4) 

1(5) 

1(6) 

1(7) 

1(8) 

1(9) 

1(10) 

1(11) 

0.8 

0.8 

0.8 

0.8 

0.8 

0.8 

0.8 

0.8 

0.8 

0.8 

0.8 

(27) 

are: 


characterizing  the  relative  flow  split  is  treated  as  object  function. 
Where  v\  is  a  function  of  l(i )  and  L(i) 

(28) 


1(1)  =  1(11),  1(2)  =  1(10), . . .  1(5)  =  1(7) 
1(1)  =  1(10),  1(2)  =  1(9),  . .  .1(5)  =  1(6) 
/(l)>/(2)>/(3)>/(4)>/(5)>/(6) 


(30) 


For  an  11 -channel  plate,  Fean  be  expressed  as: 


0.8  <  1(f),  1(f)  <  1.5 


F  = 


(29) 


Here,  F  represents  the  non-uniformity  of  the  flow  distribution.  If  F 
equals  to  zero,  all  channels  have  the  same  flow  velocity.  The  goal 
of  this  section  is  to  minimize  F  by  changing  the  variables  l{i)  and 

m 


4.3.  Constraints 

As  mentioned  above,  in  order  to  obtain  a  uniformly  distributed 
flow  field,  we  should  increase  the  cross-sectional  area  of  the  chan¬ 
nels  in  the  center  of  the  plate,  and  decrease  the  cross-sectional  area 
of  the  channels  near  the  inlet  and  outlet.  For  engineering  purpose, 
bipolar  plates  should  be  symmetric,  and  we  hope  the  width  of  chan¬ 
nels  and  ribs  are  in  the  range  of  0.8-1.5mm.  So,  the  constraints 


4.4.  Results 

The  commercial  software  MATLAB  is  used  to  optimize  the  flow 
field,  and  the  optimized  results  are  listed  in  Tables  1-3. 

Table  1  shows  the  widths  of  the  channels.  After  optimization, 
the  smallest  channel  has  a  width  of  0.8  mm,  and  the  largest  channel 
has  a  width  of  1.0983  mm.  Table  2  shows  the  length  of  the  header 
between  channels.  Table  3  shows  the  mean  velocities  of  the  chan¬ 
nels  from  analytical  model,  which  indicates  that  the  velocities  at 
different  channels  are  nearly  equal  and  a  uniformly  distributed  flow 
field  is  obtained  according  to  analytical  results. 

The  objective  parameter  F  for  the  plate  before  optimization  and 
after  optimization  is  0.0580  and  0.0065,  respectively.  Flow  distribu¬ 
tion  of  the  optimized  plate  becomes  more  uniform  than  the  original 
design.  Fig.  12  shows  the  magnitude  of  the  velocity  comparison  of 
the  contours  results  obtained  from  CFD  simulation  before  and  after 
optimization. 


Table  3 

Mean  velocities  of  the  channels  from  analytical  model  (units:  ms-1). 


*1) 

v(2) 

"(3) 

v(4) 

"(5) 

"(6) 

v(7) 

*8) 

"(9) 

"(10) 

"(11) 

0.2834 

0.2836 

0.284 

0.2842 

0.2842 

0.2842 

0.2842 

0.2842 

0.2844 

0.2848 

0.2852 

0.3" 
0.25- b 


0.2 


0.3- 
0.25- 
0.2  - 1 
0.15- 


0.15- 


0.1  - 


0.05 


Fig.  12.  Comparison  of  the  contour  results  from  CFD  simulation  (units:  ms-1 ).  (A)  Before  optimization  and  (B)  after  optimization. 
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5.  Conclusions 

The  flow  distribution  is  one  of  the  important  determinants  of 
fuel-cell  efficiency.  Unfortunately,  parallel-channel  plates  often  suf¬ 
fer  from  severe  flow  maldistribution.  In  the  present  study,  a  high 
efficacy  analytical  model,  in  the  form  of  matrix  equation,  is  pre¬ 
sented  to  calculate  the  flow  distribution  of  the  parallel  channels  of 
Z-type  configurations  and  U-type  configurations,  respectively. 

Furthermore,  the  model  is  validated  by  comparison  with  the 
results  obtained  from  CFD  simulations  and  it  has  good  agreement 
with  CFD  results. 

According  to  the  simulation  results,  it  is  found  that  a  uniform 
flow  distribution  in  the  parallel  channels  can  also  be  obtained  by 
changing  the  cross-sectional  area  of  the  channels.  Therefore,  the 
analytical  model  is  used  to  optimize  the  geometric  dimensions  of 
the  plate,  and  a  more  uniformly  distributed  flow  field  is  obtained. 
It  could  be  very  helpful  in  the  concept  design  of  flow  channel  for 
bipolar  plates. 
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